library(tidyverse)
library(modelr)
library(zoo)
library(patchwork)
library(mar)
mar <- connect_mar()
default_theme <- theme_set(theme_bw())
gogn <- read_delim("teljaragogn.csv") |>
na.omit() |>
#filter(!(ar < 1992 & vatnsfall=="Blanda")) |>
mutate(ln = log(fjoldi)) |>
group_by(vatnsfall) |>
mutate( m5 = lag(rollmean(ln, 5,fill=NA, align = "right"))) |>
ungroup()
gogn |>
ggplot(aes(ar, fjoldi, col=vatnsfall)) +
geom_point() +
geom_path() +
facet_wrap(~vatnsfall, scales = "free") +
theme(axis.text.x=element_text(angle=45, hjust=1))
gogn |>
ggplot(aes(ar, fjoldi, col=vatnsfall)) +
geom_point() +
geom_path() +
facet_wrap(~vatnsfall, scales = "free") +
theme(axis.text.x=element_text(angle=45, hjust=1)) +
scale_y_log10()
gogn |>
group_by(vatnsfall) |>
summarise(start=min(ar), stop=max(ar), fjoldi_ara=n(), medal_fjoldi=mean(fjoldi), log_m_trans=exp(mean(ln)), log_medal= mean(ln), stadalfravik= sd(fjoldi), log_stadal=sd(ln)) |>
arrange(medal_fjoldi) |>
kableExtra::kbl(digits = 2) |>
kableExtra::kable_styling("striped", full_width = F)
| vatnsfall | start | stop | fjoldi_ara | medal_fjoldi | log_m_trans | log_medal | stadalfravik | log_stadal |
|---|---|---|---|---|---|---|---|---|
| Krossá | 1998 | 2023 | 25 | 226.92 | 186.54 | 5.23 | 141.18 | 0.66 |
| Vesturdalsá | 1994 | 2023 | 25 | 229.72 | 208.91 | 5.34 | 95.63 | 0.47 |
| Búðardalsá | 2016 | 2024 | 9 | 292.33 | 251.35 | 5.53 | 163.67 | 0.60 |
| LangáSveðja | 2000 | 2023 | 24 | 962.71 | 861.83 | 6.76 | 441.17 | 0.49 |
| Blanda | 1982 | 2023 | 42 | 1518.60 | 1174.30 | 7.07 | 1059.22 | 0.78 |
| NorðuráGlanni | 2002 | 2024 | 23 | 2188.57 | 1969.47 | 7.59 | 1035.97 | 0.48 |
| Elliðaár | 1933 | 2023 | 91 | 2378.84 | 2088.49 | 7.64 | 1202.84 | 0.53 |
gogn |>
group_by(vatnsfall) |>
summarise(start=min(ar), stop=max(ar), fjoldi_ara=n(), medal_fjoldi=mean(fjoldi), log_m_trans=exp(mean(ln)), log_medal= mean(ln), stadalfravik= sd(fjoldi), log_stadal=sd(ln)) |>
ggplot(aes(log_medal, log_stadal))+
geom_point() +
geom_smooth(method = "lm")
mNorm <- lm(fjoldi~vatnsfall, data=gogn) #Meðaltöl í fjölda
mLNorm <- lm(ln~vatnsfall, data=gogn) #Meðaltöl á log
mLag <- lm(ln~vatnsfall + m5, data=gogn)
mLag2 <- lm(ln~vatnsfall + vatnsfall*m5, data=gogn)
gogn |>
spread_residuals(mNorm,mLNorm, mLag, mLag2) ->
gogn
plot(mNorm)
shapiro.test(gogn$mNorm)
##
## Shapiro-Wilk normality test
##
## data: gogn$mNorm
## W = 0.91091, p-value = 9.62e-11
plot(mLNorm)
shapiro.test(gogn$mLNorm)
##
## Shapiro-Wilk normality test
##
## data: gogn$mLNorm
## W = 0.98565, p-value = 0.01663
mLNorm <- lm(ln~vatnsfall, data=na.omit(gogn))
anova(mLNorm, mLag)
anova(mLag, mLag2)
for (v in unique(gogn$vatnsfall)){
acf(gogn$mLNorm[gogn$vatnsfall==v], main=v)
pacf(gogn$mLNorm[gogn$vatnsfall==v], main=v)
}
for (v in unique(gogn$vatnsfall)){
acf(na.omit(gogn$mLag[gogn$vatnsfall==v]), main=v)
pacf(na.omit(gogn$mLag[gogn$vatnsfall==v]), main=v)
}
summary(mLNorm)$sigma
## [1] 0.5679721
summary(mLag)$sigma
## [1] 0.5038929
listi <- c("Andakílsá"=12,"Blanda"=158,"Elliðaár"=411,"Norðurá"=1861, "Krossá"=1462,
"Langá"=1568,"Miðfjarðará"=1733, "Hítará"=956,"Sog"=2271,"Stóra-Laxá"=2339,
"Ásgarðslækur"=284,"Kúðafljót"=1485,"Héraðsvötn"=955,
"LaxáLeir" = 1589,"LaxáAðal"=1591,"Þverá"=2870,"LaxáKjós"=1590, "Álftá"=69,
"Fnjóská"=537,"Grímsá"=749,"Haffjarðará"=835,
"Hafralónsá"=851,"Haukadalsa"=883,"Hofsa"=996,"Laugardalsa"=1579,
"Svalbarðsa"=2424, "Brynjudalsá"=257,"FáskrúðDölum"=465,"Hrútafjarðará" =1106,
"HölknáÞist"=1216, "LaxáDölum"=1600,"Leirvogsá"=1615, "SeláVopn"=2080,
"StaðaráSt"=2282)
veidibok_veidibok(mar) |>
left_join(tbl_mar(mar,"ops$johannes.VEIDIBOK_TO_WATER"), by=c("vatnsfall_id"="id")) |>
left_join({fv_vatn(mar) |>
select(vatn_nr, adalvatnsfall_nr, vatnsfall_til_nr)}, by=c("water_no"="vatn_nr")) |>
left_join({fv_vatn(mar) |>
select(vatn_nafn, vatn_nr) |>
rename(adalvatnsfall = vatn_nafn)}, by= c("adalvatnsfall_nr"="vatn_nr")) |>
filter(fisktegund_id == 91, ar <2023, ar > 1973) |>
filter(adalvatnsfall_nr %in% listi |
water_no %in% listi) |>
filter(!(water_no==1568 & ar %in% 1989:1990),
!(adalvatnsfall_nr==996 & ar==2008),
vatnsfall!="Blöndulón") |>
#replace_na(list(sleppt_kodi=0)) |>
group_by(adalvatnsfall_nr) |>
collect(n=Inf) |>
mutate(vatnsfall= str_flatten(sort(unique(vatnsfall)), collapse = " ")) |>
ungroup() |>
count(vatnsfall, adalvatnsfall_nr, ar) |>
mutate(ln=log(n)) |>
arrange(adalvatnsfall_nr, ar) |>
group_by(adalvatnsfall_nr) |>
collect(n=Inf) |>
mutate( mInf = mean(ln,na.rm = T),
m5 = lag(rollmean(ln, 5,fill=NA, align = "right")),
m10 = lag(rollmean(ln, 10,fill=NA, align = "right"))) |>
ungroup() ->
test
test |>
ggplot(aes(ar,n, col=vatnsfall)) +
geom_path() +
#geom_point() +
#geom_path(aes(y=m5), linetype = 2) +
#geom_path(aes(y=m10), linetype = 3) +
facet_wrap(~vatnsfall, scales="free_y") +
theme(legend.position = "none")
test |>
ggplot(aes(ar,n, col=vatnsfall)) +
geom_path() +
#geom_point() +
#geom_path(aes(y=m5), linetype = 2) +
#geom_path(aes(y=m10), linetype = 3) +
facet_wrap(~vatnsfall, scales="free_y") +
theme(legend.position = "none") +
scale_y_log10()
test |>
group_by(vatnsfall) |>
summarise(start=min(ar), stop=max(ar), fjoldi_ara=n(), medal_fjoldi=mean(n), log_m_trans=exp(mean(ln)), log_medal= mean(ln), stadalfravik= sd(n), log_stadal=sd(ln)) |>
arrange(medal_fjoldi) |>
kableExtra::kbl(digits = 2) |>
kableExtra::kable_styling("striped", full_width = F)
| vatnsfall | start | stop | fjoldi_ara | medal_fjoldi | log_m_trans | log_medal | stadalfravik | log_stadal |
|---|---|---|---|---|---|---|---|---|
| Staðará Steingrímsf | 1974 | 2022 | 42 | 79.52 | 67.99 | 4.22 | 44.12 | 0.59 |
| Hölkná | 1974 | 2022 | 48 | 95.38 | 81.28 | 4.40 | 47.58 | 0.63 |
| Geirlandsá Grenlækur Hellisá Holtsá Hörgsá Hörgárdal Hörgsá Síðu Jónskvísl og Sýrlækur Kúðafljót og hliðarár Skaftá Tungufljót Tungulækur | 1974 | 2022 | 49 | 105.57 | 79.91 | 4.38 | 93.98 | 0.74 |
| Krossá | 1974 | 2022 | 49 | 117.43 | 96.62 | 4.57 | 73.69 | 0.65 |
| Brynjudalsá | 1975 | 2022 | 45 | 149.24 | 114.78 | 4.74 | 108.85 | 0.81 |
| Andakílsá | 1975 | 2022 | 43 | 213.98 | 168.78 | 5.13 | 177.14 | 0.65 |
| Fáskrúð | 1974 | 2022 | 46 | 225.67 | 206.44 | 5.33 | 100.18 | 0.43 |
| Svalbarðsá Svalbarðsá Þistilsfirði | 1974 | 2022 | 48 | 254.40 | 208.90 | 5.34 | 148.33 | 0.70 |
| Héraðsvötn Hofsá Vesturdal Húseyjarkvísl Miklavatn Norðurá Skagafirði Sæmundará | 1974 | 2022 | 49 | 260.96 | 208.03 | 5.34 | 170.40 | 0.71 |
| Hafralónsá Kverká | 1974 | 2022 | 47 | 280.28 | 244.80 | 5.50 | 120.35 | 0.62 |
| Bakkaá Fnjóskadal Fnjóská | 1974 | 2022 | 48 | 285.75 | 234.56 | 5.46 | 191.25 | 0.64 |
| Álftá á Mýrum Veitá | 1974 | 2022 | 49 | 286.39 | 264.62 | 5.58 | 120.30 | 0.40 |
| Laugardalsá Laugardalsvatn | 1974 | 2022 | 48 | 296.10 | 255.06 | 5.54 | 162.85 | 0.56 |
| Hrútafjarðará og Síká | 1976 | 2022 | 46 | 343.22 | 310.79 | 5.74 | 160.61 | 0.45 |
| Leirvogsá | 1974 | 2022 | 49 | 461.37 | 414.10 | 6.03 | 214.18 | 0.49 |
| Hítará | 1974 | 2022 | 49 | 479.78 | 417.89 | 6.04 | 266.10 | 0.53 |
| Haukadalsá Haukadalsvatn | 1974 | 2022 | 49 | 676.18 | 631.94 | 6.45 | 247.06 | 0.38 |
| Ásgarðslækur Búrfellslækur Sog Stóra-Laxá | 1974 | 2022 | 48 | 844.92 | 775.01 | 6.65 | 405.08 | 0.41 |
| Haffjarðará Hlíðarvatn Mýrum | 1974 | 2022 | 49 | 926.22 | 849.43 | 6.74 | 423.62 | 0.41 |
| Eyrarvatn, Þórisst. Geitabergs Laxá í Leirársveit Selós Svínadal Þórisstaðavatn Þverá Svínadal | 1974 | 2022 | 47 | 965.79 | 907.32 | 6.81 | 341.87 | 0.36 |
| Hofsá í Vopnafirði Sunnudalsá | 1974 | 2022 | 48 | 1004.00 | 844.16 | 6.74 | 539.68 | 0.65 |
| Laxá í Dölum | 1974 | 2022 | 47 | 1014.00 | 898.18 | 6.80 | 481.41 | 0.52 |
| Elliðaár Elliðavatn | 1974 | 2022 | 49 | 1055.63 | 967.19 | 6.87 | 430.47 | 0.44 |
| Selá í Vopnafirði | 1974 | 2022 | 48 | 1152.52 | 979.02 | 6.89 | 574.35 | 0.66 |
| Grímsá Borgarfirði Tunguá Lundarreykjadal | 1974 | 2022 | 49 | 1261.00 | 1178.09 | 7.07 | 445.98 | 0.39 |
| Bugða Laxá í Kjós Meðalfellsvatn | 1974 | 2022 | 49 | 1334.08 | 1238.45 | 7.12 | 540.76 | 0.40 |
| Blanda Svartá | 1974 | 2022 | 49 | 1423.14 | 1197.31 | 7.09 | 923.57 | 0.60 |
| Langá | 1974 | 2022 | 47 | 1540.55 | 1425.55 | 7.26 | 610.30 | 0.40 |
| Norðurá Borg | 1974 | 2022 | 49 | 1667.88 | 1550.12 | 7.35 | 643.39 | 0.40 |
| Austurá Miðfjarðará Núpsá Vesturá | 1974 | 2022 | 49 | 1730.14 | 1453.36 | 7.28 | 1131.13 | 0.59 |
| Kráká Laxá í Aðaldal Mýrarkvísl Reykjadalsá og Eyvindarlækur | 1974 | 2022 | 49 | 1786.92 | 1547.79 | 7.34 | 998.05 | 0.55 |
| Kjarrá Litla-Þverá Þverá Þverá Borgarfirði eystri | 1974 | 2022 | 48 | 1927.25 | 1805.30 | 7.50 | 735.15 | 0.36 |
test |>
group_by(vatnsfall) |>
summarise(start=min(ar), stop=max(ar), fjoldi_ara=n(), medal_fjoldi=mean(n), log_m_trans=exp(mean(ln)), log_medal= mean(ln), stadalfravik= sd(n), log_stadal=sd(ln)) |>
ggplot(aes(log_medal, log_stadal))+
geom_point() +
geom_smooth(method = "lm",se=FALSE) +
geom_smooth(col="black", se=FALSE)
Sleppum því að tjekka á normaldreifingu.
#Líkön
mL <- lm(ln~vatnsfall, (test))
mLm5 <- lm(ln~vatnsfall + m5, data = (test))
mLm10 <- lm(ln~vatnsfall + m10, data = test)
mLNA <- lm(ln~vatnsfall, na.omit(test))
mLm5NA <- lm(ln~vatnsfall + m5, data = na.omit(test))
mLm5m10 <- lm(ln~vatnsfall + m5+m10, data = test)
mL0m5 <- lm(ln~m5, data = na.omit(test))
mL0m5m10 <- lm(ln~m5 + m10, data = na.omit(test))
mL0m10 <- lm(ln~m10, data = na.omit(test))
mLm5NAspec <- lm(ln~vatnsfall + vatnsfall*m5, data = na.omit(test))
test |>
spread_residuals(mL, mLm5,mLm10,mLm5m10) ->
test
plot(mL)
shapiro.test(test$mL)
##
## Shapiro-Wilk normality test
##
## data: test$mL
## W = 0.9852, p-value = 2.059e-11
Feitir halar, veiðihlutfall og léleg skráning?
drop1(mLm5m10, test="F")
drop1(mLm5NAspec, test = "F")
5 ára meðaltal með mest áhrif
summary(mLm5m10)$sigma
## [1] 0.4780038
summary(mLm5)$sigma
## [1] 0.4991747
summary(mL)$sigma
## [1] 0.5475801
for (v in unique(test$vatnsfall)){
#acf(na.omit(test$m2[test$vatnsfall==v]), main=v)
pacf(na.omit(test$mLm5[test$vatnsfall==v]), main=v)
}